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^~~^ ' Abstract 

> ■ 

—^ ' The multidimensional relativistic hydrodynamical code GENESIS has 

\ jr^ I been used to obtain first results of 3D simulations of relativistic jets. We 

^+ , have studied the influence of a slight perturbation of the injection velocity 

f^ ' field on the morphodynamics of otherwise axisymmetric relativistic jets. 

G\ ■ 1 Introduction 

O .■ Astrophysical jets are continuous channels of plasma produced by some active 

' I galactic nuclei that are currently observed in radio frequences. The relativistic 

^ ■ nature of the plasma has been inferred from (esentially) two observational evi- 

"^ I dences: (i) the existence of superluminal motions of some radio components and, 

C^ ■ (ii) the high flux variability (even smaller than one day for some sources). Since 

several years the dynamical and morphological properties of axisymmetric rela- 
tivistic jets are investigated by means of relativistic hydrodynamic simulations 
rN I (e.g., [^, 0, [l5| , P^ ). In addition, relativistic MHD simulations have been per- 

formed in 2D (gg I) and 3D (|l|, 0). In their 3D simulations, |l| and [0, 
have studied mildly relativistic jets (Lorentz factor, W — 4.56) propagating both 
along and obliquely to an ambient magnetic field. 

In this work we report on high-resolution 3D simulations of relativistic jets 
with the largest beam fiow Lorentz factor performed up to now (7.09), the largest 
resolution (8 cells per beam radius), and covering the longest time evolution (75 
normalized time units; a normalized time unit is defined as the time needed for 
the jet to cross a unit length). These facts together with the high performance of 
our hydrodynamic code allowed us to study the morphology and dynamics of 3D 
relativistic jets for the first time. 

The calculations have been performed with the high-resolution 3D relativistic 
hydrodynamics code GENESIS IT] , which is an upgraded version of the code devel- 
oped by p5|. GENESIS integrates the 3D relativistic hydrodynamic equations in 
conservation form in Cartesian coordinates including an additional conservation 



equation for the beam-to-external density fraction to distinguish between beam 
and external medium fluids. The code is based on a method of lines which first 
discretizes the spatial part of the relativistic Euler equations and solves the fluxes 
using the Marquina's flux formula B . Then the semidiscrete system of ordinary 
differential equations is solved using a third order Runge-Kutta algorithm [EOl. 
High spatial accuracy is achieved by means of a PPM third order interpolation 
[pi. The computations were performed on a Cartesian domain (X,Y,Z) of size 
15i?f, X 15Rb X 75Rb (120 x 120 x 600 computational cells), where Rb is the beam 
radius. The jet is injected at z = in the direction of the positive z-axis through 
a circular nozzle defined by x^ + y^ < R^. Beam material is injected with a beam 
mass fraction / = 1, and the computational domain is initially filled with an 
external medium (/ = 0). 

We have considered a 3D model corresponding to model C2 of ||l5[, which is 
characterized by a beam-to-external proper rest-mass density ratio rj = 0.01, a 
beam Mach number Mb = 6.0, and a beam flow speed Vb = 0.99c (c is the speed 
of light) or a beam Lorentz factor Wb « 7.09. An ideal gas equation of state with 
an adiabatic exponent 7 = 5/3 describes both the jet matter and the ambient gas. 
The beam is assumed to be in pressure equilibrium with the ambient medium. 

The evolution of the jet was simulated up to T Ri 150i?;,/c, when the head of 
the jet is about to leave the grid. The scaled final time T « 4.6 lO'' (i?;,/100pc) yr 
is about two orders of magnitude smaller than the estimated ages of powerful jets. 
Hence, our simulations cannot describe the long term evolution of these sources. 

Non-axisymmetry was imposed by means of a helical velocity perturbation at 
the nozzle given by 



< = Cvb cos f j , vl ^ Cwfc sin f j , v^ = Vb\J\ - C^ 



(1) 



where C is the ratio of the toroidal to total velocity and t the perturbation period 
{i.e.,T = T/n^ n being the number of cycles completed during the whole simu- 
lation). This velocity field causes a differential rotation of the beam. The per- 
turbation is chosen such that it does not change the velocity modulus, (i.e., mass, 
momentum and energy fluxes of the beam are preserved). 

2 Morphodynamics of 3D relativistic jets 

Here we consider two models: A, which has a 1% perturbation in helical velocity 
(C = 0.01) and n = 50 and B, with ( = 0.05 and n = 15. Figure|l| shows various 
quantities of the model A in the plane y = at the end of the simulation. Two 
values of the beam mass fraction are marked by white contour levels. The beam 
structure is dominated by the imposed helical pattern with amplitudes of « Q.2Rb 
and « 1.2Rb for A and B, respectively. 

The overall jet's morphology is characterized by the presence of a highly tur- 
bulent, subsonic cocoon. The pressure distribution outside the beam is nearly 
homogeneous giving rise to a symmetric bow shock (Fig. np) in model A. Model 
B shows a very inhomogeneous pressure distribution in the cocoon. As in the 



classical case Q, the relativistic 3D simulation shows less ordered structures in 
the cocoon than the axisymmetric models. As seen from the beam mass fraction 
levels, the cocoon remains quite thin (~ 2i?h) in A and widens (^ 4i?f,) in B. 

The flow field outside the beam shows that the high velocity backflow is restric- 
ted to a small region in the vicinity of the hot spot (Fig.ll|e), the largest backflow 
velocities (~ 0.5c) being signiflcantly smaller than in 2D models. The flow annulus 
with high Lorentz factor found in axisymmetric simulations is also present, but it 
is reduced to a thin layer around the beam and possesses sub-relativistic speeds 
(~ 0.25c) in model A and mildly relativistic {^ 0.7) in B. The size of the back- 
flow velocities in the cocoon do not support relativistic beaming in case of small 
perturbations but such possibility is open in larger ones. 

Within the beam the perturbation pattern is superimposed to the conical 
shocks at about 26 and 50 Rb- The beam of A does not exhibit the strong per- 
turbations (deflection, twisting, flattening or even fllamentation) found by other 
authors Q for 3D classical hydrodynamic jets; g for 3D classical MHD jets). 
This can be taken as a sign of stability, although it can be argued that our simu- 
lation is not evolved far enough. For nl5p01, the beam is about to be disrupted 
at the end of our simulation. Obviously, the beam cross section and the internal 
conical shock structure are correlated (Figure |l|). 

The helical pattern propagates along the jet at nearly the beam speed which 
could yield to superluminal components when viewed at appropriate angles. Be- 
sides this superluminal pattern, the presence of emitting fluid elements moving at 
different velocities and orientations could lead to local variations of the apparent 
superluminal motion within the jet. This is shown in Fig.nlf, where we have com- 
puted the mean (along each line of sight, and for a viewing angle of 40 degrees) 
local apparent speed. The distribution of apparent motions is inhomogeneous and 
resembles that of the observed individual features within knots in M87 B. 

The jet can be traced continuously up to the hot spot which propagates as 
a strong shock through the ambient medium. Beam material impinges on the 
hot spot at high Lorentz factors in A case, but the beam Lorentz factor strongly 
decreases for B. We could not identify a terminal Mach disk in the flow. Wc flnd 
flow speeds near (and in) the hot spot much larger than those inferred from the 
one dimensional estimate. This fact was already noticed for 2D models by [|2| and 
suggested by them as a plausible explanation for an excess in hot spot beaming. 

We find a layer of high specific internal energy (typically more than a tenfold 
larger than that of the gas in the beam core, see Fig.|l|d) surrounding the beam like 
in previous axisymmetric models [y . The region filled by the shear layer is defined 
by 0.2 < / < 0.95. It is mainly composed of forward moving beam material at a 
speed smaller than the beam speed (Fig. |l]e) . The intermediate speed of the layer 
material is due to shear in the beam/cocoon interface, which is also responsible 
for its high specific internal energy. The shear layer broadens with distance from 
0.2Ri) near the nozzle to I.IR), (in A) or 2.0i?(, (in B) near the head of the jet. The 
diffusion of vorticity caused by numerical viscosity is responsible for the formation 
of the boundary layer. Although being caused by numerical effects (not by the 
physical mechanism of turbulent shear) the properties of PPM-based difference 
schemes are such that they can mimic turbulent flow to a certain degree |ia| . The 



existence of such a boundary layer has been invoked by several authors ml, 14 to 
interpret a number of observational trends in FRI radio sources. Such a layer will 
produce a top-bottom asymmetry due to the light aberration [EJ, and additionally, 
it can be used to explain the rails in polarization found by 0] . Other authors plj l 
have found evidence for these boundary layers in FRIIs (3C353) radio sources. 

The jet's propagation proceeds in two distinct phases. First it propagates 
according to a linear of ID phase, and then the behavior depends on the strength 
of the perturbation: it accelerates to a propagation speed which is « 20% larger 
than the corresponding ID estimate in model A or it deccelerates up to 0.37c. The 



second result partially agrees with the one obtained by [16| g/7| . 

The axial component of the momentum of the beam particles (integrated across 
the beam) along the axis decreases by more than a 30% within the first 60 Rb- Ne- 
glecting pressure and viscous effects, and assuming stationarity the axial momen- 
tum should be conserved, and hence the beam flow is decelerating. The momentum 
loss goes along with the growth of the boundary layer. 

In model A, although the beam material decelerates, its terminal Lorentz factor 
is still large enough to produce a fast jet propagation. On the other hand, in 3D, 
the beam is prone to strong perturbations which can affect the jet's head structure. 
In particular, a simple structure like a terminal Mach shock will probably not 
survive when significant 3D effects develop. It will be substituted by more complex 
structures in that case, e.g., by a Mach shock which is no longer normal to the 
beam flow and which wobbles around the instantaneous flow direction. Another 
possibility is the generation of oblique shocks near the jet head due to off-axis 
oscillations of the beam. Both possibilities will cause a less efficient deceleration 
of the beam flow at least during some epochs. At longer time scales the growth of 
3D perturbations will cause the beam to spread its momentum over a much larger 
area than that it had initially, which will efficiently reduce the jet advance speed. 
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Figure 1: Rest- mass density, pressure, flow Lorentz factor, specific internal energy and 
backflow velocity distributions (from top to bottom) of the model discussed in the text 
in the plane y = at the end of the simulation. White contour levels appearing in each 
frame correspond to values of / equal to 0.95 (inner contour; representative of the beam) 
and 0.05 (representative of the cocoon/shocked external medium interface). The panel 
under (e) displays the isosurface oi f — 0.95. Panel (f): Mean local apparent speed 
observed at an angle of 40 degrees. Arrows show the projected direction and magnitude 
of the apparent motion the contours corresponding to values of 1.0 c, 1.6 c, and 2.2 c, 
respectively. Averages have been computed along each line of sight using the emission 
coefficient as a weight. 
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